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ABSTRACT 

According to the latest evidence, the Universe is entering an era of exponential expan- 
sion, where gravitationally bound structures will get disconnected from each other, 
forming isolated 'island universes'. In this scenario, we present a theoretical criterion 
to determine the boundaries of gravitationally bound structures and a physically mo- 
tivated definition of superclusters as the largest bound structures in the Universe. We 
use the spherical collapse model self-consistently to obtain an analytical condition for 
the mean density enclosed by the last bound shell of the structure (2.36 times the 
critical density in the present Universe, assumed to be flat, with 30 per cent matter 
and 70 per cent cosmological constant, in agreement with the previous, numerical re- 
sult of Chiueh and He). iV-body simulations extended to the future show that this 
criterion, applied at the present cosmological epoch, defines a sphere that encloses 
~ 99.7 per cent of the particles that will remain bound to the structure at least until 
the scale parameter of the Universe is 100 times its present value. On the other hand, 
(28± 13) per cent of the enclosed particles are in fact not bound, so the enclosed mass 
overestimates the bound mass, in contrast with the previous, less rigorous criterion 
of, e. g., Busha and collaborators, which gave a more precise mass estimate. We also 
verify that the spherical collapse model estimate for the radial infall velocity of a shell 
enclosing a given mean density gives an accurate prediction for the velocity profile of 
infalling particles, down to very near the centre of the virialized core. 

Key words: galaxies: clusters: general - cosmology: theory - (cosmology:) large-scale 
structure of Universe 



1 INTRODUCTION 

The evidence for an accelerated expansion of the Universe, 
initially based on the observations of distant supernovae 
jGarnavich et alJ Il998l : iPerlmutter et alJ Il999l) and later 
strengthened by precise m easurements of cosm ic microwave 
background fluctuations ('Spergel et alj |2003^ , has estab- 
lished a new cosmological paradigm based on the presence 
of a 'dark energy' component. In the new scenario, the Uni- 
verse has recently made a (smooth) transition from a matter- 
dominated, decelerating stage to a dark-energy dominated, 
accelerating stage. 

In the simplest models, consistent with the observations 
so far, the dark energy behaves like Einstein's cosmologi- 
cal constant, providing an always present, constant, positive 
energy density (and a negative pressure of the same magni- 
tude). As long as the matter density was substantially larger 
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than the dark energy density, it dominated the evolution of 
the Universe, decelerating the expansion and driving the for- 
mation of structures by gravitational instability. When the 
average matter density fell below that of the dark energy, 
the latter started accelerating the expansion, and the for- 
mation of structure slowed down, as the gravitational forces 
between matter elements decreased due to their increas- 
ing separation. In this stage, structures much denser than 
the dark energy are not affected by the latter and remain 
bound, while they separate from each other at an accelerat- 
ing rate, which does not allow them to join in larger struc- 
tures. Thus, at the present cosmological time, when the ac- 
celeration of the expansion has recently started, the largest 
bound structures are just forming. In the future evolution 
of the Universe, their individual, internal properties (such 
as physical size and density) will not change substantially, 
but they will grow increasingly isolated , forming 'island uni- 
verses' (e. g., |A dams fc Laughlin 1997; Ch iueh fc h3|2002| : 
iNaeamine fc Loeb..2003. : .Busha et aL.2003il . 
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In the present Universe, these structures have not 
yet fully formed, virialized, and clearly separated from 
each other, making it difficult to identify them unambigu- 
ously. Superclusters, the largest structures identifiable in the 
present Universe, have g enerally been defined by more or less 
arbitra ry criteria (e. g., Quint ana. Carrasco fc Reiseneegeil 
[2000 ; Ei' nasto et al.l200lll2003allbHProust et all2005^ . Here. 



we address the question of how to decide whether such a 
structure will remain gravitationally bound in the future 
evolution of the Universe, forming an island universe, and 
propose to use this criterion as a physical definition of su- 
perclusters. 

We study the well-known spherical collapse model in 
the presence of a cosmological constant, with the aim of ob- 
taining a useful method to study structure behavior in our 
Universe. The spherical collapse model considers a spheri- 
cally symmetric mass distribution, where spherical shells ex- 
pand or collapse into the centre of the structure in a purely 
radial motion and without crossing each other. Specifically, 
we study the density that needs to be enclosed by a spherical 
shell at a given cosmological epoch in order to stay gravita- 
tionally attached to the overdense region within in a distant 
future, dominat ed by a cosmologic a l cons tant. 

Previously, iLokas fc Hoffmsinl i2002h gave an approx- 
imate criterion for this density, assuming that the shell 
in question initially expands with the Hubble flow (with- 
out retardation from the enclosed over-density) and present 
density evolu t ion cu rves considering different cosmologies. 
IChiueh fc Hel ||2Qq2)) on the other hand, solved the spheri- 
cal collapse equations numerically, both with a cosmolog- 
ical constant and with a more general form of dark en- 
ergy (with a constant 'equation-of-state parameter' Wq), ob- 
taining a self-consistent, theoretical criterion for the mean 
density enclosed in the last gravitationally bound or 'criti- 
cal' sh ell. Both iNagaminc & Locb (2003) and Busha ot lil] 
i2003h studied the future evolution of structures numerically, 
contrasting the extension of bound structures in a distant 
future (at scale parameter a = 166 and a = 100, re spec- 



tively) with the criterion of [Lokas fc Hoffman' f2002) ap- 
plied at the present time {a — 1). While Naeamine fc Loeb 
l|200,3r focused on the evolution of specific structures in our 
local Universe (the Local Group, the Virgo cluster, and other 
nearby structures'). iBusha et alJ |2003) followed the internal 
density and velocity structures of generic, bound objects as 
they evolve, a subject taken up a gain by the same authors 
more recently ijBusha et al.ll200,'^^ . 

We stress that, among the previous work cited above, 
only Chiuch & He (2002) used the spherical collapse model 
self-consistently in order to obtain the density criterion for 
the last bound shell, whereas all the other papers rely on the 
incorrect assumption that the shell expands with the Hubble 
flow at the initial epoch (taken to be the present time, both 
by INagaminc fc Loc b 2003 and Busha ct al. 2003). On the 
other hand, IChiueh fc Hel (|2002ll did not perform numeri- 
cal simulations of the future evolution of structure in order 
to test the accuracy of their result in a non-ideal situation 
and in order to study the details of the evolution of bound 
objects, as the other authors did. In this work, we attempt 
to bring together and extend the best parts of the previ- 
ous work, by a self-consistent application of the spherical 
collapse model (supplemented by a new, analytic equation 
for the mean density of a marginally bound sphere) and a 



comparison of its predictions to A''-body simulations of the 
future Universe. 

In section |5| we review the spherical collapse model, 
deriving an analytical solution for the spherical collapse 
equations . Our solution , whic h agrees with the numerical 
resuh of IChiueh fc Hel l|20o3), relates the critical shell's 
overdensity with the value of the dimensionless parameter 
^A{t) = A/[3H{t)^] (where A is the cosmological constant 
and H{t) is the Hubble parameter), so we are able to eval- 
uate the criterion at any time t of the evolution of the Uni- 
verse. We also obtain numerical solutions for the velocities 
of non-critical shells in the present Universe {Qa ~ 0.7). 

In sections |3| through]^ we compare the theoretical re- 
sults with simulated data, studying the quality of the crite- 
rion and its applicability to r eal-world observatio ns. As in 
iNagamin e fc Loebl (120031) and lBusha et alJ J2003h . iV-body 
simulations were run until the very late future (a = 100), 
in order to reproduce the flnal configuration of bound struc- 
tures. We find that the ri gorous density crite rion proposed 
by us (in agreement with lChiueh fc Hg|2002^ gives a good 
upper bound to the size of the bound structure, as it en- 
closes ~ 99.7 per cent of the bound particles. On the other 
hand, it also encloses a substantial number of unbound par- 
ticles, and therefore does not p erform a s well as the crite- 
rion o flLokas fc Hoffma nI (l2002h (used bv lNaeamine fc Loebl 
I2OO3I and lBushaet al.ll2003l) for the purpose of estimating 
the bound mass. We also show that the spherical collapse 
model is quite accurate in predicting (as a function of the 
enclosed mean density) the radial velocity of the stream of 
particles falling into a bound structure for the first time, 
down to the very centre of the structure. 



2 SPHERICAL COLLAPSE MODEL 

A straightforward approach to the evolution of structure 
is obtained by considering a spherical distribution where all 
layers expand or contract with only radial motions and with- 
out crossing each other. The later ensures that the enclosed 
mass in every shell stays constant throughout its whole evo- 
lution, being the only parameter other than initial condi- 
tions to determine its behavior. The model is more accurate 
during the expansion phase, but becomes unrealistic after 
contraction begins, because the instability in angular mo- 
mentum during this phase causes non-radial motions and 
finally virialization. In our analysis, we are going to leave 
aside this warning and check the accuracy of the model in 
the contraction phase using comparison with simulations. 

Our analysis w ill be based on Newtonian Mechanics, 
which, according to iLemaitrd (Il93ll) . is a limiting approx- 
imation to general relativity, valid no matter what is hap- 
pening in distant parts of the Universe. The Newtonian ap- 
proximation is accurate in a region small compared to the 
Hubble length c/H and large compared to the Schw arzschild 
radii of any collapsed objects. For more details, see iPeeble d 
Jl980l) . 

Under the assumption that the total energy will be con- 
served during the shell's expansion and later contraction, the 
evolution of a spherical shell enclosing a spherically sym- 
metric mass d istribution is given by the energy equation 
JPeeblejll980l) : 
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^ 1 /dr\2 GM A 2 



(1) 



where r is the shell's radius, M is the total mass enclosed by 
it, A is the cosmological constant and E is the total energy 
per unit mass of the shell. 

This equation can be simplified by introducing the fol- 
lowing dimensionless variables: 



^ / A 



t = 



3GM 



1/3 



(f) 



t. 



Therefore, the equation may be written 



2 \dt 



where 



E = E 



-1/3 



(2) 
(3) 

(4) 

(5) 



Fixing M and considering that there is no shell crossingj 
this equation describes the time evolution of a single shell. E 
is the dimensionless energy that describes every shell of the 
distribution, merging with the background's mean density 
when E = 0, shell that will expand with the Hubble fiow. 



2.1 The critical shell and turn-around radius 

We are looking for a critical shell that will stay at the limit 
between expanding for ever or re-collapsing into the struc- 
ture. To find the critical energy for such a shell, we define a 
potential energy to be maximized as 



V 



(6) 



The maximum of this potential occurs at r* = 1, so _£* = 
V {r") = — I is the maximum possible energy for a shell to 
remain attached to the structure. The maximum radius for 
a critical shell is 



(^) 



1/3 



so we can reinterpret the normalized radius r as 



(7) 



(8) 



2.2 Connection with the background model 

Assuming a fiat universe with cosmological constant, the 
age of the Universe (time since the Big Bang, written in 
our dimensionless variable s) may be rela ted to the vacuum 
energy density parameter JPeeblejllQSol) . 

nA = ^=tanh^{3t/2), (9) 

showing that S7a increases monotonically in time, and can 
therefore be used as a surrogate time variable, in terms of 
which we will describe the evolution of bound or critical 
(marginally bound) spherical shells. 

Integrating equation Q from the beginning of time (r = 
0) till the current radius of a given shell (r = ro), we get 




Figure 1. The normalized radius res of a critical shell as a func- 
tion of Ha- The dotted lines highlight the present universe, in 
which Ca = 0.7 and res = 0.84. 



to 



Vrdr 



y/r^ + 2Er + 2 



(10) 



In the particular case of a critical shell {E — E* = —^), the 
denominator of equation llOt can be factorized, yielding 



to = 



'•0 



Vrd/r 



,0 il-7)Vr + 2 
The integral can be done analytically, with the result 



to 



2^3 



In 



l+2ro + y/Sro (rp + 2) 
l + 2n,- y/3n-> {ro + 2) 



(11) 



(12) 



In 



1 + ro + y^ro (r„ + 2) 



Notice that ro can be chosen at any 'current' time, so 
equation l|12|l can be generalized to any time by simply drop- 
ping subscript 0. From now on, we will use the subscript 'cs' 
to indicate that we are referring to a critical shell. 

Replacing eq. H12|l in eq.Q, and introducing a new vari- 
able 



X(rcs) = 



1 -I- 2rcs + y/Svcs (res + 2) 
1 + 2Fes - y^Srcs (res + 2) 

< (^1 + Fes + V^cs (res + 2) 



(13) 



we can write the relation between Qa (and therefore cosmo- 
logical time) and res as 



i^A (res) = 



X(rc; 



X(rcs) + 1 



(14) 



As expected (see Fig. 0, res grows with Qa as 57a grows 
with time, and will converge to its maximum radius when 
t — » oo. 

Inversely, evaluating equation (I14II at our preferred cos- 
mology (f2ni — 0.3, SIa = 0.7), we obtain that current critical 
shells should have res ~ 0.84, which means that their present 
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radius is 84 per cent of the maximum radius they will reach 
as t ^ <x. 



2.3 Conditions for a critical shell 

For practical applications, it is convenient to express the 
critical condition as the minimum enclosed mean density 
needed by a shell to stay bound to the central gravitational 
attractor. The critical density of the Universe (in the flat 
model, also the average total density, including matter and 
vacuum energy) is 



Pc = 



8^' 



(15) 



and the average mass density enclosed by a given shell is 



Pm = 



3A/ 



Defining the mass density parameter for the shell. 



(16) 



(17) 




Figure 2. Density parameter for a critical shell, £7s (see equa- 
tion 1171 ). as a function of Qa and dimensionless radius Tcs- Its 
maximum value, f2cs = 2.37, occurs at Ha = 0.72 (res = 0.85). 



the condition for the shell to be bound is 



2.36, 



(18) 



where we have evaluated equation 11411 at the present value 
of Q,A = 0.7. This represents the present density contrast in- 
side the last shell that will eventually stop its growth at the 
end of times. This result was first obtained bv .Chiueh fc If a 
i2002h . using numerical methods, but now we confirm this 
result with an analytical solution. In the same way, evalu- 
ating at Q.A = 1 (t oo) we obtain the asymptotic critical 
density contrast condition 57cs,oo = 2. 

Fig-Hlshows the density parameter fJcs as a function of 
Q.h and res using equation 1141 or its inverse. It is curious 
that the curve shows a single maximum at S7a = 0.72, very 
close to the measured value of S7a today. Considering that 
f^A can be taken as a function of time, this means that we 
are living in the era when the ratio between mass density 
inside a critical shell and the critical density of the Universe 
is very nearly at its maximum. 

For observational purposes, it is more interesting to 
know the ratio between the mass density enclosed by the 
critical shell, p^, and that of the background, p°^. Using 
equation 1171 1. we obtain that 



^e 



Pm Pc(l-f^A) I-^Ja' 

which, evaluated for Q,a ~ 0.7, yields 



^ = 7.i 



(19) 



In the notation of lBusha et alJ i2003h . considering only 
the excess density with respect to the background (Mobj ~ 
M — Mb, where Mb is the amount of mass contained in 
an equivalent sphere with background density) and writ- 
ing the present value of the Hubble parameter as Ho ~ 
70 hrokm s^^ Mpc~^, we rewrite the criterion as 



ro 



IMpc 



(20) 



This is a less restri ctive condition than the one proposed by 
iBusha et all J2003l) . which is 



M, 



obj 



1012 Mp 



3 hro 



1 Mpc 



(21) 



or, equivalently, Qs ^ 5.56 (also equivalent to the overden- 
sity criterion given b y Naeaminc & Locb 2003, based on the 
formalism of iLokas fc Hoffman .2002.) . Thi s result was ex- 
pected, since the criterion of I'Busha et all ll2003f) was ob- 
tained under the assumption that the test particles, placed 
on the critical shell today, move with the Hubble flow. As 
this is not really true, they free that parameter (velocity of 
particles on the critical shell today) and do an empirical test 
to observe the real behavior of particles evolved till a = 100. 
They obtain a corrected result which is closer to their initial 
result than to ours, always setting a higher constraint to the 
critical shell density. A probable reason for this difference is 
that our result is a completely theoretical approach based 
on the spherical collapse model (only radial motions), but 
in real life objects usually obtain angular momentum caused 
by tidal forces, which will tend to detach objects from the 
structure, strengthening the binding condition. This discus- 
sion will be resumed later, in relation to the results of our 
simulations. 



2.4 Velocities of shells 

In order to work with shell velocities, it is convenient to refer 
the radial velocities to the Hubble flow, so we introduce the 
parameter A, deflned as in iBusha et all l|20o3), 

/ 1 dr^" 

For a critical shell, we can use equation Q to write this in 



(22) 



terms of 17a as 



-Acs {^A 



+ 



r2 (J7a) r3 (Oa) 



(23) 
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Figure 3. The dimensionless radius r of a shell as a function of 
its dimensionless energy i? jit a cosmological time characterized 
by = 0.7. A shell with E = —2.1 corresponds to r = 0, so it 
is now collapsing. Highlighted by dotted vertical and horizontal 
lines are the critical shell, with Ec = —1.5 and Vc = 0.84, and the 
shell at turn-around, with i?ta = —1.64 and rta = 0.73. 

Evaluating at Qa ~ 0.7, we obtain that the present velocity 
parameter for a critical shell is Acb = 8.63 x 10~^, show- 
ing that the shell has been slowed down substantially with 
respect to the Hubble flow. 

To study the velocity profile of a mass overdensity, it is 
convenient to relate the energy of an arbitrary shell to its 
normalized radius, characterized by the present value of Qa- 
Unfortunately, the energy integral has no analytical solution 
for E ^ E* , so this relation can only be obtained numer- 
ically. Fixing the value of to 0.7, in order to represent 
the present universe, we obtain the current normalized time 
to — 0.81. Then, we numerically integrate equation llOII for 
every possible value of E from zero till some value of r that 
satisfies the time constraint. Now using equation 117L we ar- 
rive at the more useful numerical function E{^ls; ^a), which 
is best shown in Fig. |3| 

Expressing the velocity parameter A in terms of an ar- 
bitrary normalized energy E{Q,s',^a), for a fixed value Qa, 
we obtain 

A = nA + Qs + E(ns; Qa) (2J7a!:!s . (24) 

In Fig. 21 we observe that A{Qs) is zero when Qs ~ Qs,ta., 
the density parameter corresponding to the turn-around ra- 
dius. When solving for the radial velocity, it is important 
to notice that shells with fis > ns,ta contract, while shells 
with fis < f2s,ta expand. For a universe with Qa = 0.7, 
fis.ta — 3.66. 



3 NUMERICAL SIMULATIONS 

We simulate one cosmological model, assuming a stan- 
dard flat Lambda Cold Dark Matter (ACDM) Universe. 
The current cosmological parameters in the simulation are 
f^m.o = 0.3, flA,o = 0.7 and h = 0.7, where the Hub- 
ble parameter Ho — lOOft km s~^Mpc~^. The normaliza- 
tion of the power spectrum is cts = 1. The box has a side 




Figure 4. Radial velocity parameter A (see eq. 1241 for a shell as 
a function of the enclosed density Qb (eq. 1171 in a universe with 
f^A = 0.7. For a critical shell, Ucs = 2.36 and A(Oa) = 0.086. For 
a shell at its turn-around point, A = and f2s = 3.66. 

length of 100/i~^Mpc and contains 128'^ dark matter par- 
ticles of mass TTlDM = 3.97 X lO^^/i'^M©. The simulation 
was evolved from a = 0.02 (redshift z = 49) to a = 100. 
The Plummer-equivalent gravitational softening was set to 
epi = 15 h^^ kpc (physical units) from a — 1/3 to a = 100, 
while it was taken to be fixed in comoving units at higher 
redshift. 

The run was performed wi th the massive parallel tree 
7V-bo dy/SPH code GADGET iSpringel. Yoshida fc Whitd 
I2OOII) . This is a TREESPH code where the dark mat- 
ter particles are evolved using a tree-code, while the col- 
lisional gas is followed using the SPH approach^. Here, 
we used the new improved GADGET2, kindly provided by 
Volker Springel dSpringell lioOS) . which is more memory- 
efficient, offers better time-stepping for collisionless dy- 
namics, and is substantially faster than the original ver- 
sion. The initial conditions were establ ished by the code of 
Ivan de Weveaert fc BertschingeJ il996l) . 

We took snapshots at the present time (a — 1) and 
in the far future (a — 100), assuming that in late epochs 
the structure evolution will decrease significantly, so no ma- 
jor changes will be seen from then on. The identification of 
structures was done using a friend-of-friends code to identify 
the initial candidates and later extracting the structures to 
produce a reduced catalog for our study. 



4 APPLICATION OF BINDING CRITERIA TO 
SIMULATION DATA 

We identified 22 mass concentrations with masses greater 
than 10^* h^gMQ, from which we selected 11 that were 
away from the box boundaries in order to do our analysis^. 
The same identification was done in the future frame (a = 

^ The gas particles were not considered in our simulations. 

We could have used the other 11 objects as well, because the 
box had periodic boundary conditions, but we did not, just to 
simplify the analysis. 
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Figure 5. Enclosed mass as a function of critical radius for 11 
bound objects identified in the simulation at a = 1. The dotted 
line represents the constraint of 2.36 times the critical density, 
used to define the bound region. 



100), obtaining the final state for tlie structures identified 
at a = 1. We were able to follow particles from the present 
till the late future frame, so we could exactly determine the 
fate of every particle. 

In order to apply the density criterion for gravitational 
binding, we first identified the densest core of the struc- 
ture. For this, we chose a centre and found the radius where 
the mean inner density was 300 times the critical density 
of the Universe (chosen to identify a dense, virialized core 
that would be clearly detectable in observations). We then 
re-centred the sphere to its centre of mass. We repeated 
this procedure until the centre of mass matched the geo- 
metric centre of the sphere. Once the centre was fixed, we 
calculated the density parameter fig of concentric spheres 
with increasing radius until the condition fig < fics = 2.36 
stopped being satisfied. We called the identified radius Tc- 
(The radii and masses of these structures are plotted in Fig. 
|SJ|. The same identification was done at a = 100, with the 
critical condition f2a < fics.oo = 2. I n order to contrast our 
results, we applied the criterion of iBusha et alJ (|200^ at 
a = 1, which is given by Q.^ < ilca,B = 5.56. We called the 
identified radius vb- 

A critical part of the analysis is to have an acceptable 
criterion to determine whether a particle is bound to the 
overdense structure or not. As illustrated in Fig.|S| the the- 
oretical criterion at a = ICQ gives a very intuitive result, 
placed at the end of the virialized region, almost exactly 
where the lower envelope of the radial velocities crosses zero. 
In other words, Slcs.oo coincides with the last radius where 
objects with negative radial velocity can be found. At that 
point, only a few particles will be able to escape, so we can 
say that the criterion at a = 100 is adequate to determine 
the limits of bound structures in the distant future universe. 
A very similar figure, as well as a detailed discussion of how 
the velo city distribut i on evo lves to this state, was recently 
given bv lBusha et"all (|200^. 

Finally, we classified particles in four categories: parti- 
cles that fall inside the theoretical critical shell at a = 1 and 




log{r/kpc) 

Figure 6. Radial velocities as a function of radius (in physical 
units) for a mass overdensity at a = 100. The vertical dashed line 
shows the radius where Qb = f2cs oo = 2. 



a = 100; particles that fall inside the critical shell at a = 1 
but do not at a = 100; particles that fall outside the criti- 
cal shell at a = 1 but inside at a = 100 and particles that 
fall outside in both cases. This categorization lets us visual- 
ize the quality of our estimation, clearly separating particles 
according to their real fate, and permitting us to calculate 
statistical indicators useful to produce the desired error es- 
timations. As expected, not every object presented a very 
'spherical' distribution at a = 1, and this affected how well 
data were fitted by the spherical model''. For this reason, we 
based our qualitative analysis on the 'best' objects, but we 
kept all of them when doing statistics. 



5 COMPARISON BETWEEN THEORY AND 
SIMULATED DATA 

The spherical model predicts a purely radial motion of par- 
ticles towards the centre. This is clearly not the case in the 
real world, where objects are affected by multiple and com- 
plicated tidal forces produced by other objects in their sur- 
roundings all along their evolution. In fact, objects present 
a velocity dispersion that is greatest in the virialized cores 
of galaxy clusters. The presence of other overdensities sur- 
rounding the main attractor alters the motion, pushing par- 
ticles away from their radial trajectories. For this reason, 
we expect the spherical model to give a lower bound on the 
radial velocities at a given radius (with outward-pointing ve- 
locities taken as positive) and an upper limit for the radius 
of critical shells. 

To test the performance of the criterion, we first selected 
all the particles that were bound at a = 100 and 'marked' 
them, so they could be recognized at a = 1, as shown for the 
most massive structure in Fig. Q Later, we selected all the 
particles that satisfied the criterion at a = 1 and produced 



^ We even found some objects that were currently undergoing 
mergers, so they showed strong evolution between a = 1 and 



a ■■ 



100. 
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Figure 7. Spatial distributio n of particles in object no. 1 (see Table [Tl with colours depending on our criterion (left-hand panel), or 
on that of lBusha et al.l i20n3h (right-hand panel). Blue; they satisfy the criterion at a = 1, but not at a = 100, Green: they satisfy the 
criterion at a = 1 and a = 100, Red: they do not satisfy the criterion at a = 1, but they do at a = 100. Particles that do not satisfy the 
criterion at either a = 1 or a = 100 are not shown. 



four statistical indicators^: A, the fraction of all particles 
selected by the criterion as bound at a = 1 that were not 
selected at a = 100; B, the fraction of all particles selected 
at a = 1 that were also selected at a = 100; C, the ratio of 
the number of particles not selected at a = 1 but selected at 
a — 100, to the total number of selected particles at a = 1 
and D, the ratio of the mass still bound at a = 100 to 
the mass selected at a = 1. Results were averaged over all 
objects to obtain a single indicator. Results for our criter ion 
are shown in Tableland for that of lBusha et alj i2003l) in 
Tabled 

We observe that the purely theoretical criterion is good 
to determine the exterior limit of the object since only very 
few particles (0.26 per cent of the particles that satisfied the 
criterion at a = 1) that are outside the criterion today will 
fall inside in the late future. In contrast, a significant number 
of particles (28 per cent of the particles under the criterion) 
were predicted to belong, but finally escaped. Finally, a large 
number of particles (72 per cent of the particles under the 
criterion) were correctly predicted to belong to the overden- 
sity^. Based on this result, we can assert that the theoretical 
criterion is adequate to give an upper bound to the struc- 
ture extension, but overestimates its mass, which turns out 
to be about 72 per cent of t he predicted mass. 

With the criterion of iBusha et al.l (|200^, a greater 
number of particles (13 per cent of the particles under the 
criterion at a = 1) fell outside the criterion today, but ended 
up inside the object. In contrast, a smaller number of parti- 



^ Only two of these are independent, as A -|- _B = 1 (or 100 per 
cent) and B + C = D. 

^ Mean values from the 11 objects studied. 



cles (10 per cent of the particles that satisfied the criterion) 
were incorrectly predicted to belong, while a very large num- 
ber of particles (90 per cent of the particles under the crite- 
rion) were correctly recognized. Thus, although it is based on 
an inconsistent application of the spherical collapse model, 
it happens to give a better estimate for the final object mass, 
underestimating it by only ~ 3 per cent. 

It is also important to add that objects 6 and 7 had 
very complex cores, being separated into two main concen- 
trations. Of these, we selected the most massive as the centre 
of the spherical analysis. We kept our selection criterion for 
the centre, under the assumption that it should be easier to 
detect the centre of the most massive concentration when 
dealing with real observations. 



5.1 Radial velocity predictions 

An important prediction of the spherical collapse model is 
the radial velocity of shells falling toward the gravitational 
attractor. The velocity information is contained in the ve- 
locity parameter A, which depends only on the mass energy 
density inside the shell, Og. This can be calculated numeri- 
cally for every shell of the studied objects, yielding the de- 
sired radial velocity. 

In Fig.|H| we plot the radial velocity against radius, to- 
gether with the theoretical approximation using the spheri- 
cal collapse model with and without cosmological constant 
(see Reise ncgger et aL 200Q,, for details about the approx- 
imation using the spherical collapse model without a cos- 
mological constant). For this analysis, we chose object no. 9 
because it presented a clear stream of infalling particles fiow- 
ing at high speed to the object's core, cleanly separated from 
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Table 1. Quantitative results after application of our theoretical criterion to 11 objects from the simulation. Columns: (1) index of the 
object; (2) A, particles selected by the criterion as bound at a = 1 but not at a = 100, as fraction of all particles selected at a = 1; (3) 
B, particles selected at a = 1 and a = 100, as fraction of all particles selected at a = 1 (note A + B = 100 per cent); (4) C, ratio of 
the number of particles not selected at a = 1 but selected at a = 100, to the total number of selected particles at a = 1; (5) M(a = 1), 
total mass selected inside the critical radius at a = 1; (6) M(a = 100) total mass selected inside the critical radius at a = 100 and (7) 
D = M(a = 100)/M(a = 1). 



Object 


A 


B 


C 


M{a = 1) 


M(a = 100) 


D 


no. 


[per cent] 


[per cent] 


[per cent] 


[10" h-^ Mq] 


[10" h-^ Mq] 


[per cent] 


1 


16.8 


83.2 


0.20 


7.63 


6.36 


83.4 


2 


18.3 


81.7 


0.22 


6.46 


5.29 


81.9 


3 


18.7 


81.3 


0.50 


5.76 


4.71 


81.8 


4 


37.9 


62.1 


0.04 


6.66 


4.14 


62.1 


5 


25.4 


74.6 


0.62 


5.41 


4.07 


75.2 


6 


51.6 


48.4 


0.21 


6.64 


3.23 


48.6 


7 


50.5 


49.5 


0.06 


5.96 


2.96 


49.6 


8 


26.7 


73.3 


0.07 


3.34 


2.45 


73.3 


9 


17.2 


82.8 


0.67 


2.88 


2.40 


83.4 


10 


17.6 


82.4 


0.09 


2.45 


2.02 


82.4 


11 


29.6 


70.4 


0.22 


2.07 


1.46 


70.6 


Mean 


28.2 


71.8 


0.26 






72.0 


Std. dev. 


13.0 


13.0 


0.23 






13.1 



Table 2. Quantitative results after application of the criterion of lBusha et al.l 1200.^) to 11 objects from the simulation. The variables 
are defined as in Table [Tl 



Object 


A 


B 


C 


M(a = 1) 


A/(a = 100) 


D 


no. 


[per cent] 


[per cent] 


[per cent] 


[10" fe-l Mq] 


[10" /i-i Mq] 


[per cent] 


1 


7.4 


92.6 


9.3 


6.24 


6.36 


101.9 


2 


8.7 


91.3 


6.0 


5.44 


5.29 


97.2 


3 


7.1 


92.9 


27.7 


3.91 


4.71 


120.5 


4 


10.8 


89.2 


2.7 


4.51 


4.14 


91.8 


5 


8.4 


91.6 


17.2 


3.74 


4.07 


108.8 


6 


10.2 


89.8 


42.8 


2.43 


3.23 


132.6 


7 


15.8 


84.2 


5.9 


3.28 


2.96 


90.1 


8 


16.3 


83.7 


7.2 


2.69 


2.45 


90.8 


9 


3.7 


96.3 


18.8 


2.09 


2.40 


115.2 


10 


5.9 


94.1 


4.8 


2.04 


2.02 


98.9 


11 


14.9 


85.1 


4.3 


1.64 


1.46 


89.4 


Mean 


9.9 


90.1 


13.3 






103.4 


Std. dev. 


4.2 


4.2 


12.4 






14.3 



the virialized particles and from particles currently outflow- 
ing after a first pass through the centre, forming an empty 
space reaching as close as 2Mpc from the core. Other ob- 
jects presented similar overall characteristics, but they did 
not show such a clear separation between infalling particles 
and virialized ones, probably as a result of earlier virializa- 
tion or the presence of substructure. We observe that the 
spherical model with cosmological constant is accurate to 
determine the mean radial velocity of the infalling part of 
the cluster particles, correcting the underestimate of radial 
velocities at large radius seen in the spherical model without 
cosmological constant. 

An important observation is that the theoretical veloc- 
ity profile follows the infalling particles deep into the core 
of the structure, where virialization effects are very impor- 
tant. This result, which was confirmed on every object we 
studied, tells us that the spherical collapse model produces 
robust predictions of the negative velocity envelope profile 



even in highly virialized cores. The decrease on the predicted 
velocity very near to the centre is probably due to the poor 
resolution of the simulation and to extreme virialization ef- 
fects. 



5.2 Perturbations from spherical collapse 

In the previous analyses, we observe that there is a consid- 
erable number of particles that escape from the structure, 
contradicting the theoretical prediction. There are two main 
candidates to be responsible for this: one is the appearance 
of tangential velocities during the contraction, which is also 
responsible for virialization, the other is the influence of ex- 
ternal structures which can act as gravitational attractor on 
external shells. A simple test is to plot tangential velocities 
as a function of radius. We would expect that if there is a 
clear relation between angular momentum and failure of the 
criterion, we would find that particles that contradict the 
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Figure 8. Radial velocities as a function of radius for a mass overdensity at a = 1 for object no. 9. Left-hand panel: colours as a function 
of tangential velocities. Red dots indicate more radial velocity and bluer ones less radial velocity. Right-hand panel: colours depending 
on final fate of particles: green particles remain bound to the structure, blue particles do not. The vertical dashed line shows the radius 
where Qs < ^cs = 2.36, while the vertical dotted line indicates the critical radius according to lBusha et alJ 1200,'^) . The solid curve is the 
theoretical approximation considering a flat universe with a cosmological constant. The dashed curve is the theoretical approximation 
considering a universe without cosmological constant and f2j\/ = 0.3. 
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Figure 9. Tangential velocities as a function of radius for object 
no. 9 at a = 1. Colors depend on the final fate of particles, as 
given in Fig. 171 



criterion would have greater tangential velocities than par- 
ticles that verify the criterion. This was done for object no. 
9 and results are shown in Fig. 1^1 

As seen in Fig. |U1 and as we observed in all other ob- 
jects, there is no clear relation between tangential velocity 
and probability of escaping. This result is solid enough to 
say that angular momentum is not the main reason for the 
failure of the criterion, so the hypothesis of external attrac- 
tors seems more acceptable. Just as a follow-up observation, 
we can see in Fig. |7] that the particles incorrectly identified 
as bound by our criterion (blue) are commonly related to 
denser regions outside and close to the critical shell, while 
the particles incorrectly identified as unbound (red) are close 
to mass concentrations inside the critical shell that finally 
fall to the centre. We conclude that the failure of the method 



is mainly caused by external overdensities or perturbations 
from the spherical distribution of the studied object. 



6 CONCLUSIONS 

We have presented a complete discussion of the model of 
spherical, gravitational collapse in a flat Universe with a cos- 
mological constant, applied to estimate the size and mass of 
bound structures in the Universe. Within this model, we de- 
rive an exact, analytical equation for the minimum enclosed 
mass density required for a shell to remain bound until a 
very distant future. For the present cosmological parameters 
(fim = 0.3, Qa = 0.7), this minimum density is 2.36 times 
th e critical dens i ty of the Universe, as found numerically 
bv lChiueh fcH^ (|2003)- This suggests a both physical and 
practical criterion for the limits of a bound structure, such 
as the previously fairly ill-defined superclusters of galaxies, 
as the shell enclosing precisely this density. 

The application of the model to simulated data gave 
encouraging results, demonstrating first its great ability to 
find the limits of structure in the late future, and, second, 
giving reasonably good results in determining the limits of 
bound structures at the present epoch. On average, 72 per 
cent of the mass enclosed by the estimated radius is really 
bound to the structure, while the mass that, although bound 
to the structure, is not enclosed by the radius is only a 0.3 
per cent. For the less rigorous criterion of Busha ct al. ( 20oJ) 
and other authors, these numbers are 90 and 13 per cent, 
giving a substantially better estimate for the object's final 
mass. Thus, the sphere defined by our criterion is an outer 
envelope enclosing all the particles bound to the structure 
(and quite a few more), while that of .Busha et al. (200.J) 
encloses as much mass as will remain bound to the distant 
future (leaving about as many bound particles outside as 
unbound ones inside). 

The spherical collapse model also defines a radial veloc- 



10 R. Diinner et al. 



ity profile that we will use in our next paper (Diinner et al., van de Weygaert R., Bertschinger E., 1996, MNRAS, 281, 

in preparation) to find the shape of these structures in red- 84 
shift space, in order to make them identifiable in redshift sur- 
veys. This profile was found in the simulations to fit the ob- 
served velocity profile of infalling particles well down to deep 
inside the virialized radius. Thus, even though the spherical 
collapse model is intrinsically unstable for contracting shells, 
it still gives a reliable performance in a broad set of radii. 
Moreover, we observed that the greatest perturbations from 
the theoretical model were produced by gravitational per- 
turbations by external structures or by substructures falling 
into the external shells of the structure. 
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